Solvent-pumped evaporation concentration on paper in linear and radial geometries

Solvent-pumped evaporation-driven concentration of an initial distribution of solutes on a porous substrate is considered in one and two dimensions. Approximate analytic solutions to the isotropic advection–dispersion equations are first found for a Gaussian kernel and an infinite domain, following the smoothed particle approximation. Analytic solutions for more general initial distributions are then found as sums of Gaussians, and comparison is made with numerical solutions. In each case, initial distributions are advected toward the stagnation point and concentrated. Two-dimensional analysis is then extended to describe anisotropy in permeability and diffusion, and hydrodynamic dispersion. Radial-flow experiments are performed using filter papers and water-soluble dyes. Diffusion coefficients, temperature and humidity profiles, and the evolution of spot distributions are measured. The results confirm minor anisotropy in permeability and diffusion, limited hydrodynamic dispersion, and largely uniform evaporation. Péclet numbers over 2500 are demonstrated. Evaporation-driven concentration provides a mechanism for solute transport over long timescales. Potential applications lie in the design of paper spray microanalytical devices operating by solvent pumping rather than capillary flow.

In a previous paper, 12 we compared the performance of 1D and 2D evaporation-driven paper concentrators with infinite sources and highlighted the advantage of radial flow. Advective transport toward a single stagnation point provides an effective counter to diffusion. In devices of comparable Péclet numbers (which defines the relative significance of the two effects), 37 radial flow offers faster concentration than linear flow after a longer filling time and may provide the basis for other 2D microfluidic devices. However, one difficulty preventing application may be the problem of finding analytic solutions for ADEs with varying coefficients, which masks potential capabilities.
Analytic solutions for special cases have been found in hydrology. [38][39][40] More general initial value problems have been solved using power series 41 and Laplace transforms, 36,42 if necessary, using numerical inversion. 43 As a last resort, ADEs can always be solved numerically. 44,45 An alternative approach is offered by smoothed particle hydrodynamics (SPH), which provides approximate solutions in terms of assemblies of kernel functions. Originally developed in astrophysics, 46 SPH has also been successful in fluid flow problems involving porous and anisotropic media. [47][48][49] Here, we use an SPH approach to model solventpumped concentration of solute distributions in 1D and 2D porous systems, with the overall aim of understanding transport, concentration, and coalescence of an initial distribution of spots.
The paper is arranged as follows. In Sec. II, the ADEs for 1D and 2D geometries are reviewed. In Sec. III, SPH solutions are developed for Gaussian initial distributions in 1D domains and extended to general distributions. In Sec. IV, analysis is repeated for 2D domains. In Sec. V, the modifications needed to model anisotropy in permeability and diffusion are introduced. In Sec. VI, hydrodynamic dispersion is considered via numerical solutions. In Sec. VII, radial-flow experiments with dyes on filter papers are presented. Measurements of diffusion coefficients and temperature and humidity profiles are used to validate assumptions of isotropic diffusion and uniform evaporation. The evolution of initial distributions is then demonstrated, confirming anisotropic permeability but little hydrodynamic dispersion even with fan-assisted evaporation. Applications are considered in Sec. VIII, and conclusions are drawn in Sec. XI.

II. ADVECTION-DISPERSION EQUATIONS
The ADEs for evaporation-concentration on a porous substrate have been derived elsewhere. 12 Here, we give a brief outline. Figure 1 shows the geometries of (a) linear and (b) radial systems. In each case, a substrate of porosity ε, thickness d, and finite dimension (X M or R M ) is assumed to contain an initial distribution of non-volatile solute, which is carried by a volatile solvent from a source S to a concentration point C, with solvent evaporation between. In the linear case, the source is distributed over a width w, and in the radial case, it is distributed round the perimeter. We assume that the substrates are horizontal (so gravity may be neglected) and that the capillary filling phase has ended (so that substrates are fully wetted). We also initially assume isotropy and uniform evaporation; these assumptions will be reviewed later. For a 1D system, Darcy's law 50,51 implies that the local solvent velocity V satisfies Here, P is the local fluid pressure, k is the permeability, and μ is the dynamic viscosity. Assuming a uniform evaporation flux v e per unit area from both sides together, continuity implies that the solvent flow velocity satisfies Here, τ e ¼ d/v e is a characteristic time constant. Combining the above, we obtain k μ Assuming V ¼ 0 at X ¼ 0 and P ¼ P atm at X ¼ X M , where P atm is atmospheric pressure, we obtain The velocity, therefore, falls linearly toward the origin, while the pressure falls parabolically; it is this pressure drop that provides pumping.
The ADE can be derived by considering the balance of solute transport, diffusion, and accumulation. For a linear adsorption isotherm, the resulting equation can be written in terms of normalized quantities as @C Lx @t 0 ¼ 1 Pe Here, C Lx (x) is the liquid phase solute concentration, M /Dτ e is the Péclet number, and D is the diffusion coefficient, initially assumed to be both isotropic and constant. Finally, t 0 ¼ R f T/τ e is the normalized time, where T is the actual time and R f ¼ 1/(1 þ aF) is the retardation factor. Here, a is Henry's coefficient and F ¼ (1 À ε)/ε is the volumetric ratio of the stationary and mobile phases. The left-hand side of the ADE describes accumulation, and the right-hand side models diffusion and convective transport.
For a radially symmetric 2D system, it is simple to show that the velocity profile and pressure distribution modify to V ¼ ÀR/2τ e and P ¼ P atm À μ(R 2 M À R 2 )/4kτ e , so the peak velocity and pressure drop are both halved compared to a 1D system with equivalent dimensions. For a radially symmetric concentration profile C Lr (r), the ADE can then be written in terms of normalized  r @C Lr @r þ r 2 @C Lr @r þ C Lr : Here, the Péclet number is now Pe ¼ R 2 M /Dτ e . However, we are interested in the evolution of non-axisymmetric concentration profiles C Lxy (x, y) in the same radial flow field. In this case, the ADE can be written in Cartesian coordinates as @C Lxy @t 0 ¼ 1 Pe Questionable aspects above are the neglect of capillary transport, and the assumptions of uniform evaporation, isotropic permeability, and a constant, isotropic diffusion coefficient. Before proceeding, we make the following comments: Capillary transport has been considered by many authors, both without 52-54 and with [55][56][57] evaporation. Its effect is to redistribute a solute deposited on a dry substrate. Once wetting is complete, evaporation will concentrate the new distribution by solvent pumping, and this step is the focus here. The evaporation rate will be determined by kinetics; again, this aspect has been considered by many authors. 58,59 The main effect, convection, will be uniform if the paper temperature is uniform. Variations are likely to be small near the stagnation point, especially in the radial case. However, the evaporation rate will be different near the source, where heat transfer also takes place by conduction. This effect will be considered further in Sec. VII. Paper has preferred axes (the machine and cross directions). Anisotropic permeability leads to curved flow trajectories. 60,61 Anisotropy may also alter the diffusion tensor. 62,63 These effects will be considered further in Sec. V. Hydrodynamic dispersion in porous media leads to non-Fickian diffusion described by velocity-dependent axial and transverse diffusion coefficients. 49,64,65 Because of the low flow velocity, these effects are likely to be limited here, but will be considered further in Sec. VI.

III. 1D INITIAL VALUE PROBLEMS
We begin with initial value problems in one dimension. To do so, we first assume a trial solution to the ADE in the form of a moving Gaussian kernel, as Here x 0 and f x are functions of t 0 that describe the time evolution of the center and width of the Gaussian, and C L0 and f x0 are constants. Differentiation, substitution into the ADE and equation of coefficients of terms (x À x 0 ) n separately then gives Both these equations may be solved by integration. For the initial conditions x 0 (0) ¼ x 00 and f x (0) ¼ f x0 , the results are Because x 0 tends to zero, all initial distributions will evolve to a Gaussian at the origin. The solutions all tend from f x0 at t 0 ¼ 0 to unity as t 0 ! 1, and the steady state solution is The solutions above are exact for infinite space, but do not satisfy the boundary conditions for a finite concentrator. To obtain a suitable solution, we assume mirrored initial conditions. In this case, the analytic solution becomes By symmetry, the solution derivative must now be zero at x ¼ 0. If the Gaussian terms always decay rapidly enough, the modified solution can be a good approximation for a finite concentrator over the reduced domain 0 x 1. Figure 2(a) compares the analytic solution (blue), with a numerical solution of the ADE (green), for Pe ¼ 1000 and an initial Gaussian distribution with a half-width and position defined by Pe 0 ¼ f x0 Pe ¼ 250 and x 0 ¼ 0:5. The solute profile initially travels to the left, narrowing and increasing in height as it does. When it reaches the origin, it is trapped, and gradually evolves into a half-Gaussian. Figure 2(b) shows a similar comparison for the time variation of the concentration at x ¼ 0. Its value is initially zero, when the Gaussian is some distance from the origin. As it approaches the origin, it rises rapidly and then tends to a steady state. In each case, there is excellent agreement between the analytic and numerical solutions.
Because the ADE is linear, solutions for arbitrary initial distributions can be constructed as sums of Gaussians. Unfortunately, the expansion coefficients cannot be determined uniquely. However, reasonable approximations are possible if care is taken to manage the overlap of adjacent terms. Good results can be obtained using sums of identical, offset Gaussians, as with Here, the terms x i are a set of N þ 1 centers with separation Δ ¼ 1/N, and the terms a i are expansion coefficients chosen to match the initial conditions, as The factor of 2 is used to prevent double counting at x i ¼ 0 and the weighting factor K is introduced to compensate for the overlap of adjacent Gaussian terms, as A smooth approximation may be obtained by assuming a suitable value for the initial Péclet number Pe 0 ¼ f x0 Pe of each Gaussian. If Pe 0 is too small, and the centers are too widely spaced, the approximate distribution will contain ripples. However, if Pe 0 is too large, it will be difficult to model fast-varying distributions. Good results are obtained with Pe 0 ¼ 2N 2 . In this case,  and their amplitudes lie some way below the target profile because of the factor K. However, the overall sum yields a realistic approximation to the target profile. Its flat top is recreated smoothly; its edges are modeled less abruptly but become sharper if N is increased. Figure 3(b) shows the predicted concentration profile at the different times t 0 in the legend (blue); these agree almost exactly with a finite difference solution of the ADE for the same initial distribution (green).

IV. 2D INITIAL VALUE PROBLEMS
We now repeat the process for 2D initial value problems. We first assume a trial solution to the ADE for a single offset Gaussian kernel as Here x 0 , y 0 , f x , and f y are all functions of t 0 and describe the evolution of the Gaussian center and shape, and C L0 , f x0 , and f y0 are constants. For the initial conditions f These solutions are analogous to previous 1D results; but time constants have been doubled by the halved flow velocity. Once again, the solutions describe the evolution of an initial distribution to a steady state at the stagnation point, but no extra solution components are required to satisfy boundary conditions. Provided the Gaussian is not too extended, the solution is again a good approximation over the domain 0 r 1.
The left-hand equations imply that motion is always toward the origin. Similarly, the right-hand equations imply that the spot shape gradually evolves, with the limiting values f x1 ¼ f y1 ¼ 1 showing that the eventual shape is circular. For example, Fig. 4 shows concentration profiles calculated using the analytic solution at the different times t 0 in the legend, for Pe ¼ 1000, and an offset elliptical Gaussian initial distribution with Pe 0x ¼ 20Pe, Pe 0y ¼ Pe/20, x 0 ¼ 0:8 and y 0 ¼ 0. As t 0 increases, the distribution gradually tends to the origin and becomes more and more circular. We have verified that there is excellent agreement with numerical solution of the ADE for the same initial distribution. It is simple to show that the steady state solution is Comparison with Eq. (14) shows that the peak width is doubled in the 2D case.
The procedure above may easily be adapted to model the evolution of a set of Gaussian spots with different sizes, shapes, and centers, using one term for each spot. A different approach is again required for arbitrary distributions. In this case, a solution may be constructed as a sum of displaced circular Gaussians, by assuming with : (26) Here, the coordinates x ij and y ij are a grid of centers with separation Δ ¼ 1/N, and the terms a ij are coefficients chosen to match the initial condition C L0 (x, y) as A smooth approximation may again be obtained by assuming a suitable value for the initial Péclet number of each Gaussian. We have obtained good results with Pe 0 ¼ 4N 2 . Numerical summation of the Gaussians contributing to each point then shows that K % π: (28) Figure 5(a) shows line-scans in the x-direction of an expansion for a disk-shaped distribution extending up to r ¼ 0:6, with Pe ¼ 1000. Here, we choose N ¼ 50 and Pe 0 ¼ 10000. Individual Gaussians are shown in blue, their sum in green, and the target distribution in red. Once again, adjacent Gaussians overlap, and their amplitudes lie below the target profile because of the factor K. However, the sum yields a realistic approximation. Figure 5(b) shows the predicted concentration profiles at the different times t 0 in the legend (blue); these agree almost exactly with a finite difference solution of the ADE for the same initial distribution (green). The corresponding y-variations are identical. Figure 6 shows 2D plots of the concentration profile at different times t 0 , which clearly evolves from a uniform disk to a centered Gaussian peak.

V. ANISOTROPY
We now revisit the assumption of isotropy. If the Cartesian axes are aligned with the principal directions of an anisotropic substrate, solvent flow may be described by velocities V x (X) and V y (Y). Continuity implies that In the same way, Darcy's law may be taken for an anisotropic substrate as 60,61 Here, k x and k y are the permeabilities in the two directions. Combining the above we obtain Biomicrofluidics ARTICLE pubs.aip.org/aip/bmf A symmetric solution for the pressure with stagnation at R ¼ 0 and P ¼ P atm at R ¼ R M can be found as This result again implies that the pressure falls parabolically toward the origin, but that the velocity components are now different, We can write these alternatively as Here, α is a coefficient that describes the flow anisotropy, namely, Clearly, α ¼ 0 represents isotropy. Macroscopic solvent streamlines then have the form Integrating from (X 0 , Y 0 ), we then have For isotropic permeability (α ¼ 0), all trajectories are straight; for anisotropic permeability, oblique trajectories are curved. This observation allows estimation of anisotropy in permeability from spot trajectories.
For anisotropic media, we may also assume that diffusion is described by different diffusion coefficients D x and D y , 62,63 written conveniently as D x ¼ βD and D y ¼ β À1 D, where D ¼ ffiffiffiffiffiffiffiffiffiffi ffi D x D y p and β ¼ 1 represents isotropy. The ADE can again be derived by considering the balance of solute transport, diffusion, and accumulation, as A solution for a single offset elliptical Gaussian spot may again be attempted in the form of Eq. (21), which for convenience we repeat here, Here C L0 , f x0 , and f y0 are again constants, and the functions x 0 , y 0 , f x , and f y now satisfy For the initial conditions x 0 (0) ¼ x 00 , y 0 (0) ¼ y 00 , f x (0) ¼ f x0 , and f y (0) ¼ f y0 , the solutions are Biomicrofluidics ARTICLE pubs.aip.org/aip/bmf Here, the constants γ and δ are given by We now consider the implications, starting with peak position. The center of the Gaussian spot again tends to the origin, and the time variations of x 0 and y 0 are simply different exponentials. However, the trajectory slope is now Solute spot trajectories therefore follow solvent streamlines. We now consider the peak shape. In each case, the solution rises

Biomicrofluidics ARTICLE
pubs.aip.org/aip/bmf from the starting value to a limiting steady state, given by If β ¼ 1, results for f y are obtained by reversing the sign of α in results for f x . Generally, the final spot shape depends on both α and β. For example, f x1 is increased (and the corresponding Gaussian width is decreased) if β , 1 (so that D x is relatively low) or if α . 1 (so that V x is high). However, an initially circular spot will remain so despite anisotropy if γ ¼ δ, when differences in diffusion are exactly balanced by differences in advection. This requires the anisotropic coefficients to be related by

Biomicrofluidics ARTICLE
pubs.aip.org/aip/bmf Pe 0 ¼ 250, assuming that Pe ¼ 2500. The superimposed black lines show spot trajectories. In each case, the three spots are gradually transported toward the center, concentrating as they move. In Fig. 7, α ¼ 0, β ¼ 0:5, so permeability is isotropic, and all trajectories are straight. However, diffusion is anisotropic, so the spots gradually become elliptical with their long axes in the y-direction. In Fig. 8, α ¼ À0:25, β ¼ 1 so diffusion is isotropic, but permeability is anisotropic. Oblique trajectories now become curved and spots elliptical with their long axes in the x-direction. These results imply that diffusion anisotropy can be estimated from spot shapes if isotropy in permeability is already known.

VI. HYDRODYNAMIC DISPERSION
We now briefly consider the possibility of hydrodynamic dispersion, assuming anisotropic permeability but for simplicity neglecting diffusion anisotropy. In this case, the components of the Here, α L and α T are coefficients of longitudinal and transverse dispersion, and V ¼ V x i þ V y j is the velocity vector. We can also define the normalized quantities, Here, β L,T ¼ αL,T jVj Djvj and v ¼ v x i þ v y j is a velocity vector normalized to the maximum. The ADE may then be written as This equation has no obvious analytic solution but can be solved numerically. Figure 9 shows a finite-difference solution for the three-spot pattern in Fig. 8, assuming that α ¼ À0:25, β L ¼ 16, and β T ¼ 0 (anisotropic permeability and large hydrodynamic dispersion). As before, trajectories become curved. However, spots are no longer elliptical, and their long axes rotate away from principal directions toward streamlines, implying that hydrodynamic dispersion can be differentiated from diffusion anisotropy by alteration and reorientation of spot shapes.

VII. EXPERIMENTAL RESULTS
Radial concentration was performed with Brilliant Blue FCF, a water-soluble dye. Its diffusion coefficient has been measured in water as D ¼ 5:68 Â 10 À10 m 2 /s. 66 Dye was obtained at ≥97% purity (80717, Sigma Aldrich, St. Louis, USA) and dissolved in de-ionized water. Concentration was carried out on 90 mm diameter qualitative circles (1001090, Whatman, Little Chalfont, UK, thickness d ¼ 180 μm, pore size 11 μm), on which the retardation factor of Brilliant Blue is ≈1. 67 Experiments were performed in laboratory conditions . Despite this variation, results were qualitatively repeatable. Figure 10 shows the arrangement used. The paper was held horizontally between two plastic supports to leave a portion of radius R M ¼ 30 mm suspended. The lower support contained an annular reservoir linked to the evaporation section, while the upper support acted as a paper-cutting guide. Water drops were dispensed to locate the principal axes. Radial slots were scalpel-cut to allow tabs to be bent into the reservoir, which was filled with de-ionized water, so the paper was fully wetted. Dye spots were then dispensed. Evaporation could be prevented by sealing the plastic jig with transparent lids, performed under ambient conditions, or enhanced by using a bladeless fan (X12, Archuu, Shenzhen, China).
Diffusion, transport, and concentration were monitored using a Dino-Lite AM2011 Basic USB digital microscope (AnMo Electronics Corp., Taipei, Taiwan), arranged to view through the fan duct, with a white InGaN-on-sapphire LED strip-light (LuckyLight, Shenzhen, China) for illumination. Figure 11 shows the experimental jig with the microscope and fan in place. Spot shapes and trajectories were extracted by post-processing time-lapse photographs to extract transmission profiles from red-channel data and converting to spatial variations in concentration using the Beer-Lambert law. Although the light source is broadband, linearity has been established for the concentration range here. 12 The diffusion coefficient of Brilliant Blue on wetted papers was found by photographing a single dye spot at different times T, with evaporation suppressed. Best-fits between line-scan absorbance data and Gaussian profiles were then obtained in the form Values of c were then extracted, and the diffusion coefficient D was estimated by comparison with the solution of Fick's second law, 68 which implies that c ¼ ffiffiffiffiffiffiffiffiffi 4DT p . Figures 12(a) and 12(b) show example absorbance profiles in the cross and machine directions at the times T shown in minutes, compared with Gaussian fits. Figure 12(c) shows the variations of c with ffiffiffi ffi T p for the machine and cross directions, compared with straight line fits. Here, the slopes are almost the same, but an average of six data sets yielded D ¼ 3:08 + 0:11 Â 10 À10 m 2 /s (machine), and D ¼ 2:31 + 0:1 Â 10 À10 m 2 /s (cross) implying minor diffusion anisotropy. The values are around 45% of the earlier literature value. 66 This reduction is due to the tortuosity T of the paper microstructure, which is, in turn, related to the porosity ε. Many models show an effective diffusion constant of the form D 0 ¼ D/T 2 . One review 69 yields values of D 0 /D in the range of 0.25-0.4 for ε ¼ 0:5, slightly lower than our result. However, the experiments are hard to perform repeatedly. Similarly, diffusion anisotropy is expected from anisotropy in tortuosity.
For ambient evaporation, local temperature and humidity were measured together using a probe constructed from a data logger based on a capacitive humidity sensor and a thermistor (Thermo Recorder TR-72Ui, T&D, Matsumoto City, Japan). After removal from its original housing, the sensor head was mounted in a capsule with a small opening. Air was drawn slowly through the capsule by an axial flow fan, while the probe was translated across the paper at a small (1 À 2 mm) separation from the surface. The inset of Fig. 13 shows the sensor elements and sampling probe. The local partial pressure P W of water was estimated from the saturation vapor pressure P S and the percentage humidity f.
Here, θ is the dry-bulb temperature in degrees Celsius, and P ¼ 1018 mbar is the absolute pressure (measured here at a local weather station). P w was then found from humidity measurements, using the relation Figure 13 shows variations of θ and P W along a paper diameter, obtained as an average of three data sets. The former falls toward the center of the disk, consistent with evaporative cooling, while the latter rises. Near the center, conditions are approximately constant, confirming the earlier assumption of uniform evaporation. The ambient temperature was 22.6°C, so the central decrease is 1.2°C. This value is lower than expected because the thermistor and sampled air both reach thermal equilibrium before each measurement. Although efforts were made to provide insulation, the equilibrium temperature is higher than the inlet air temperature. Furthermore, the extrema measurements are not expected to reach ambient temperature because the probe diameter prevents the inlet being placed directly over the perimeter, and the probe then averages over several mm. For fan-assisted evaporation, temperature measurements were obtained using a FLIR E60 uncooled microbolometer array (FLIR Systems, Wilsonville, OR). The inset of Fig. 14 shows the camera and the main figure shows temperature profiles along the diameter in perpendicular directions, with and without fan assistance.
In each case, the variation is again bowl shaped, with a nearuniform central section. However, the temperature clearly falls further with fan assistance, and the region of constant temperature increases. The inset also shows a colorized thermal image obtained under these conditions. This contains a large, radially symmetric region of uniform low temperature, further confirming the assumption of uniform evaporation.
For ambient evaporation, Péclet numbers were small (∼200) and experiments took many (>10) hours. Large spot sizes complicated spot tracking and degraded the agreement with exponential models. Much larger (>2500) Péclet numbers were obtained with

Biomicrofluidics ARTICLE
pubs.aip.org/aip/bmf fan assistance, and agreement with theory improved significantly. Figure 15 shows a three-spot distribution at different times (in min) during a typical fan-assisted experiment. The spots travel from the periphery to the center of the disk and coalesce after ∼2.5 h. Radial spot extension indicates minor hydrodynamic dispersion, but its effects disappear near the stagnation point. Figure 16 shows the log-scale time variation of normalized X-and Y-positions for horizontal and vertical trajectories without and with fan assistance, extracted using a peak-finding routine and matched to straight lines with slopes a x and a y using a least squares non-linear algorithm. There is a good match to an exponential variation in each case, and the small difference in slope implies minor anisotropic permeability. The inset shows spot trajectories in the fan-assisted case.
Parameters were extracted as follows. The permeability ratio was estimated as k y /k x ¼ a y /a x , the time constant as τ e ¼ R f /(a x þ a y ), the rate of solvent evaporation as v e ¼ d/τ e , and the Péclet number as Pe ¼  Table I. There is some variation in permeability ratio, broadly in agreement with published data for paper. 60 The Péclet number is increased 13-fold by fan-assisted evaporation, rendering timescales realistic for laboratory work.

VIII. APPLICATIONS
Paper microfluidics is now highly developed. However, almost all applications are driven by capillary flow; evaporation is considered rarely. Exceptions include paper spray, 71 where extended solvent supply and desolvation have been adopted, 72,73 and reaction acceleration on surfaces and in droplets. [74][75][76] The previous results show that evaporation pumping can provide concentration as well as transport, especially in radial flow, and can operate over long timescales. Applications are likely to lie in similar areas, as follows:   (1) Concentration: The 2D analysis implies that single droplets will be concentrated at the center of a disk, or, by extension, at the tip of a wedge as shown in Fig. 17(a). The effect must lead to the sensitivity enhancement seen in paper spray mass spectrometry with solvent evaporation, allowing the best use of a limited volume of solute. (2) Improved paper substrates: The relatively long timescales of evaporation-driven concentration are a limiting factor. However, several techniques are available to manipulate permeability and evaporation by surface patterning. We are exploring laser cutting and have already used the radial-flow apparatus to demonstrate engineered orthotropic anisotropy in permeability. Materials of this type may be used to construct paper spray sources with more effective concentration as shown in Fig. 17(b). (3) Non-linear effects: The 1D analysis implies that single droplets will be concentrated at the center of a paper strip fed from a reservoir at either end. An analyte and a reagent may then be concentrated simultaneously, so that (for example) an individual concentration factor of 10 2 leads to an enhancement of 10 4 in the reaction rate. Arrangements may be combined with paper spray as shown in Fig. 17(c).
Optimizing such applications will require numerical solution of ADEs for more complex geometries and flow fields. However, the initial analysis and experimental results presented here can already allow anticipation of likely effects.

IX. CONCLUSIONS
Evaporation-driven concentration on porous substrates has been considered in 1D and 2D, and analytic solutions have been found to the ADEs by assuming uniform evaporation and infinite domains. Solutions are constructed from Gaussian kernels following the smoothed particle approximation, which provides an ideal basis for modeling spot trajectories. Provided the Gaussians are well confined, these remain valid over finite domains. The results show that initial distributions always form Gaussian concentration peaks, and that timescales are doubled in the radial geometry. The effect of anisotropic permeability is to modify spot trajectories, while anisotropic diffusion and hydrodynamic dispersion alter spot shapes. Anisotropy can be compensated by operating on the bias or eliminated using isotropic substrates, and hydrodynamic dispersion is unimportant near the stagnation point.
Theoretical predictions are confirmed by radial-flow experiments with water-soluble dyes on filter papers. Further work is required to model and directly measure local evaporation rates, but measurements of temperature, humidity, and spot trajectories all support an assumption of approximately uniform evaporation over most of the paper. Measurements have shown minor anisotropy in both diffusion and permeability. Timescales have been reduced using fan assistance, allowing transport from the perimeter to the center of a 60 mm diameter disk in a few hours. Péclet numbers over 2500 have been achieved, with limited hydrodynamic dispersion. It is likely that faster concentration can be achieved with more powerful airflow, and that the concentration point may be moved using directed air jets.
Radial evaporation-driven concentration may provide a simple method of reconcentrating forensic samples, preconcentrating samples to improve sensitivity in paper spray mass spectrometry or concentrating analytes before reaction with targeted reagents. The advantage over capillary pumping is higher sensitivity and continued operation long after capillary filling has completed. Similarly, the advantage over field-driven pumping is the simplicity of the substrate and lack of any requirement for high voltages. However, environmental control would be needed for fully repeatable results.
Relatively little solvent is required. The necessary volume V SE can be estimated as Ignoring anisotropy in permeability, spot trajectories can be written in terms of V SE as To reduce the argument of the exponential to a given value ÀC, we then require Thus, the exact evaporation rate is unimportant; V SE is dependent only on the substrate dimensions and the solute retardation factor and commensurate with the volume V SF ¼ πdR 2 M needed for filling. Assuming d, R M , and R f are as previously given, and C ¼ 4 (say), we then obtain V SE % 4cc. However, larger volumes would be required for solutes with lower R f . Transport times also increase as R f reduces.